function [hat_beta, SE, reject, reject0, pval, pval0] = implement_test(dy, z, shares, shifters, critical , ols)

%This script implements simplified version of the test used in simulations
%Corresponding author: Rodrigo Adao
%Date: 09/11/2024

Nobs = length(dy);

%Statistic estimation
if ols == 1
   %OLS version -- normalize covariance by variance of IV
   hat_beta = ( (z'*dy)/Nobs ) /( sum(z.^2)/Nobs );  
   hat_eta = dy - hat_beta*z;
else
   %Covariance version 
   hat_beta = (z'*dy)/Nobs ;
   hat_eta = dy - ( hat_beta/( sum(z.^2)/Nobs ) )*z;
end  
hat_eta0 = dy; 

%Variance with residuals computed under the null)
R0 = hat_eta0'*shares / Nobs; 
hatV0_beta = (R0.^2)*(shifters.^2); 
      
%Variance with  residuals computed with estimated OLS coefficient)  
R = hat_eta'*shares / Nobs; 
hatV_beta = (R.^2)*(shifters.^2);

%Standard error and tests
SE  = ( hatV_beta.^(1/2)  );
SE0 = ( hatV0_beta.^(1/2) );
if ols == 1
    SE = SE/( sum(z.^2)/Nobs );
    SE0 = SE0/( sum(z.^2)/Nobs );
end

%rejection decision and p-value
tstat  = hat_beta/SE;
tstat0 = hat_beta/SE0;

reject  = abs(tstat)  > critical;
reject0 = abs(tstat0) > critical;

pval  = 2*(1 - normcdf(abs(tstat) ,0,1));
pval0 = 2*(1 - normcdf(abs(tstat0),0,1));

end